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Abstract 

The physics of high-energy collider experiments asks for delicate com- 
parisons between theoretical predictions and experimental data. Signals 
and potential backgrounds for new physics have to be predicted at suffi- 
cient accuracy. The accuracy as well as the computational complexity of 
the calculations leading to the predictions depend on both the number of 
external particles in the process analyzed and the order of the quantum 
corrections, the number of loops, included in the calculation. We present 
some approaches to problems occurring in these calculations regarding the 
integration of phase-space and the inclusion of one-loop corrections. 

1 Introduction 

Collider experiments played and continue to play a fundamental role in particle 
physics. This is exemplified by past, current and future experiments like, e.g. 
LEP at CERN, Tevatron at Fermilab or HERA at DESY, and LHC, again at 
CERN, respectively. Current research focuses on even more precise tests of what 
is known as the Standard Model of particle physics, the validity of which needs 
the discovery of the Higgs boson at such a collider experiment. This asks for new 
colliders to reach higher energies leading to events from scattering experiments 
with higher numbers of particles involved. 

In order to prepare and analyze the outcome of the experiments, signals and 
potential backgrounds for new physics have to be predicted at sufficient accuracy. 
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Most of the calculations involved are set up within the framework of perturbation 
theory, in which the accuracy is controlled by the order parameter. In the appli- 
cation to quantum field theory, on which the Standard Model is based, the order 
parameter is connected to both the number of external particles in the process 
analyzed and the order of the quantum corrections, the number of loops, included 
in the calculation. 

The Standard Model is tested mainly by the comparison of the experimental 
data and theoretical calculations on the statistical level. This can happen through 
the analysis of cross sections, or simulation of collider experiments. Both the 
increase in the multiplicity of particles at the experiments, and the need for higher 
accuracy of theoretical calculations lead to a dramatic increase of the complexity 
of this analysis. In the following, we will encounter two particular examples in 
which this increase of complexity constitutes a challenge in the scientific process. 

2 A hierarchical phase space generator for QCD 
antennas 

The reliable description of multi-jet production at collider experiments is an im- 
portant issue in the study of the Standard Model. It requires the calculation of 
cross sections, which again requires the integration of squared scattering ampli- 
tudes over phase space. The number of dimensions of the integration space in 
combination with the desired cuts point at the Monte Carlo method as the only 
suitable candidate for this task. Since the computation of QCD scattering matrix 
elements with many particles is rather time-consuming, the integration process 
should preferably involve as few integration points as possible. The strong peak- 
ing structures exhibited in the QCD amplitudes enforce the application of the 
Monte Carlo method to be dressed with a sophisticated portion of importance 
sampling. Flat phase space generators, like RAMBO [1], will not be adequate for 
this task. 

In the last years several methods to efficiently integrate the peaking structures 
of the scattering amplitudes have emerged, and have been used in several con- 
texts [2]. For instance, PHEGAS [3] is an example where an efficient, automated, 
mapping of all possible peaking structures of a given scattering process has been 
established. The algorithm is based on the "natural" mappings dictated by the 
Feynman graphs contributing to the given process, so that the number of kine- 
matic channels used to generate the phase space is equal to the number of Feyn- 
man graphs. Using adaptive methods, like multi-channel optimization [4] and 
by throwing away channels that are negligible, we may end up with a few chan- 
nel generator exhibiting high efficiency, as is indeed the case in n(+7)-fermion 
production in e + e~ collisions. In contrast, the QCD scattering amplitudes point 
towards the opposite direction: large number of Feynman graphs which means 
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large number of kinematic channels which, moreover, contribute equally to the 
result. 

A way out off this problem may be based on the long-standing remark that 
n + 2-gluon amplitude may be described by a very compact expression when 
special helicities are assigned to the gluons, which, combined with the leading 
color approximation, results to 

ElMI 2 = sM (^c-i) Efe-^) 4 E A n+2 ( Pl ,..., Pn+2 ) , (f) 

c V Z 7 l<i<j P(2,...,n+2) 

where N c refers to the number of colors, 

A n+2 (pi,.. . ,p n +2) = [ (Pi ■ P2)(P2 ■ Ps) ■ ■ ■ (Pn+l ■ Pn+2)(Pn+2 ■ Pi) }^ , (2) 

and the sum over all permutations of the 2 nd to the (n + 2) nd argument of this 
function is taken, with the exception of those that are equivalent under reflection 
i i— > n + 4 — i [5] . 

SARGE [6] is the first known example of a phase space generator that deals with 
the momentum structures entering the above expression, namely with (2), known 
as antenna structures. The algorithm is based on the "democratic" strategy to 
generate the n body phase space, as is the case for RAMBO, and it makes use of 
the scale symmetry of the antenna to achieve the required goal. Now, we study 
the "hierarchical" strategy for phase space generation in order to efficiently map 
the momentum antenna structures. The idea is as follows. Using the standard 
two-body phase space (neglecting factors of 2n) 

d^ 2 (P;s 1 ,s 2 ;pi,p 2 ) = d 4 pi<J+(Pi- si)rf 4 j52^+(P2 - s 2 )5 4 (P - Pi - p 2 ) , (3) 
we decompose the phase space 

d<S> n (P; Pl ...,p n )= 5 4 (f>-P)nd 4 Pi MP*-^) (4) 



i=i 



as 



d$ n (P;pi . . . ,p n ) = ds n -id® 2 (Q n ;a n ,s n -i;p n ,Q n -i) 

X ds n - 2 d$ 2 (<2n-i; Cn-1, S n _ 2 ;p n -1, Qn-2) 

x ds 2 d$ 2 (Q 3 ;a 3 ,s 2 ;p 3 ,Q 2 ) 

x d$ 2 (Q 2 ;a 2 ,a 1 ;p 2 ,p 1 ) . (5) 

The task is to express the phase space in terms of the invariants pi -pj appearing in 
the antenna structure (2), so that, using a suitable mapping, we can construct a 
density that, apart from constant and soft terms, will be identical to this antenna 
structure. 
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2.1 Antenna generation 

In this section, we will present a sketch of how the antenna structures can be 
generated using the hierarchical approach. For a treatment in full detail, we refer 
to [7]. The starting point is the generation of he two-body phase space (3) in 
terms of the variables 

9i - Pi P2 ■ 92 , R v 

where 91,52 are given massless momenta. Let us introduce the notation 

c = cos(Z(5 1 ,5 2 )) , s = Vl-c 2 (7) 

and 

s = P 2 , s 1>2 = si, 2 /s . (8) 
We find that the parameterization 

Pi <- (S + Si - 8 2 )l{2y/s) , 

Pi *~ Pi-Vsai , 

P 2 i <- {{V~s-pl-V~sa2)+cp\)/s , (9) 
P\ - ±( (p?) 2 - si - (pi) 2 - {p\f )^ , 

leads to the identity 

d<S> 2 (P;s 1 ,s 2 ;pi,p 2 ) = da 1 rfa 2 n(a 1 ,a 2 )~ 1/2 e(n(ai,a 2 )) , (10) 

where 

Il(ai, a 2 ) = 4s 2 [(l - a 2 + s 2 - si)a 2 - s 2 ] 

- [(1 - 2ai - si + s 2 ) + (1 - 2a 2 - si + s 2 )c] 2 , (11) 

and is the step function. So in order to obtain a two-body phase space with 
a density which depends on the invariants ai,a 2 following some given function 
f(ai, a 2 ), one has to generate ai, a 2 following a density proportional to /(ai, a 2 ) x 
n(ai, a 2 ) _1//2 in the region where IT(ai,a 2 ) > 0, and construct the momenta as 
given above. 

The generation strategy proceeds through a sequence of two-body phase space 
generations following the decomposition (5). At each two-body generation, one 
final-state momentum pk is generated, together with the sum Qk-i of the remain- 
ing final-state momenta to be generated. This suggest to label the momenta in 
a way opposite to the order of generation, so first p n , Q n -i are generated, then 
Pn-i, Qn-2 and so on. The starting point is the center of mass frame (CMF) of the 
initial momenta 91 and g 2 with Q n = 91 + 92 being the overall momentum. The 
CMF of momentum Qk we denote by CMFfc. The pair pk,Qk-i is generated by 



4 



generating variables a[ k \ and constructing the momenta as described before. 
These variables are now equal to 

(k) Pk+l-Pk , (k) q2'Qk-l /10 x 

a\ = — and a> 2 = — - . (12) 

Pk+i ■ Wk <?2 ■ Vfc 

This happens in CMF^, so in order to obtain p k , Qk-i, the constructed momenta 
have to be boosted such that (\fQl, 0,0, 0) is transformed to Q k . Now, we do 
three observations. Firstly, we have 

Pk + i-Q k = (Ql +1 -Qt-P 2 k+1 )/2 ■ (13) 

Secondly, we have, with S fc = Yh=i 

Sk — Sfe d ( s k _i — S fc _x \ 

log , (14) 



(fifc — cr*: — Sfc-i) (sfc-i — Hfc_i) rfsfc_i \s k — a k — s k -i 
and thirdly, we can write 

A n+2 (qi,p n ,Pn-i ■ ■ ■ ,Pi,Q2) = ^y(s« - S„)(gi • Q„)(g 2 • <5n) 

x (n^-O^rw) -w~m > ( 15 ) 

where g fc (sfc_i) is given by (14), and where s n = Q 2 n , p n+1 = q 1 and Q 1 = p x . 
These observations suggest that the phase space generation 

ds n _i# n (s„_i) da^ da { 2 l) -^ n~ 1 ) /2 e(n (n) ) 

a\' a 2 ' 

ds n _ 2 g n -i{s n - 2 ) da { ^ 1] ^4 n_1) ^rry n (i»-?)6( n (n-i)) 

a l a 2 

ds 2 g 3 (s 2 ) daf ] da { 2 3) n^ /2 e(n (3) ) 

rfa f) J_ rfa ( 2 ) J_ n (2 ; /2 e(n (2) ) , (16) 

a\ a 2 ' 

will lead to a density for the momenta that is proportional to A n+2 . Three 
variables are generated in each CMF fc , namely Sfc_i, and a 2 k \ Just as the 
integration of s k _i (14), also the integration of results in a volume factor 

that depends on the corresponding variables generated in CMF^i. However, 
these factors are logarithmic functions of their arguments and exhibit a non- 
singular behavior, and we call them soft factors. The total actual density will 
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therefore be the product of n — 1 soft factors times the antenna structure under 
consideration. 

In the end, we want to generate all permutations in the momenta of (15). 
Those for which qi and q 2 each appear in two factors (none of which is q± ■ q 2 ) 
cannot be obtained by simple re-labeling. In order to obtain these, we observe 
that they can be decomposed into two antennas, namely 

A m+2 (qi,p m , ■ ■ - ,Pl,?2) x A n _ m+2 (q 2 ,p n , ■ ■ ■ ,p m +i,qi) (17) 

and each of these can be generated after the decomposition, 

d$ n (P;pi...,p n ) = ds m ds 

n—m d& 2 (Q n , ^mi Sn—mi Qmi Qn—m ) (18) 
n—m (Q n— rn'i Pm+1 >•••> Pn) ■ 

In order to combine the two sub-antennas to the required antenna structure, we 
have to take into account in the first decomposition a density that is proportional 
to 

© ( \/ $n \f yJSn—m ) 



(qi ■ Q m )(qi ■ Q n—m )(?2 • Qm){q2 ■ Qn-m) S m S n - 

2.2 Results 



(19) 



In this section, we present results obtained by SARGE and HAAG 1 , the program 
that implements the hierarchical algorithm described before. In order to be as 
general as possible, the only cut we apply is 

fe+^) 2 >so , (20) 

where i, j(i ^ j) runs from 1 to n+2 where n is the number of final-state particles. 
Unless explicitly mentioned differently, we use sq = 900 GeV 2 and the total energy 
y/s = 1000 GeV. Moreover, all particles are assumed to be massless in order to 
compare with SARGE, with which only massless particles can be treated. 

As it was mentioned in the introduction, we are interested in integrating sums 
of QCD antenna structures (2). We start by considering the simplest case, namely 
integrating the function 

S 2 [ (Pi -Ps)(P3 -Pa)(P4, -P2)(P2 -Ph) ■ ■ ■ (Pn+2-Pl) }^ (21) 

that corresponds to a given permutation of the momenta, namely (1,3,4,2,5, 
. . . , n + 2). In Table 1 we give the results for SARGE, and HAAG. In all three codes 
the same single channel, corresponding to (21), has been used in the generation. 
iV gen and iV acc are the number of generated and accepted events, and by / we 
define 

f = f 2 , (22) 



1 HAAG stands for: Hierarchical AntennA Generation. 
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where V2 is the quadratic variance and / is the estimated integral. / is clearly a 
measure of the efficiency of the generator. Moreover e, defined as 



'"max 

is the usual generation efficiency related for instance to 'unweighted' events in a 
realistic simulation. The results agree well, and exhibit the fact that the generated 
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SARGE 
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13.6 


0.23 
1.02 


6 


SARGE 
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1 x 10° 
1.2 x 10 5 
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.487 x 10" 11 


.164 x 10~ 12 
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1141 
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0.21 
1.48 



Table 1: Results for the single-channel integration/generation. 



densities of the generators the hierarchical type are much closer to the integrand. 
The same picture is reproduced for an arbitrary permutation. 

For a realistic QCD calculation, the integrated function may be approximated 
by a sum over permutations. Therefore, an efficient generator has to include all 
possible channels, where each channel corresponds to a given permutation of 
the momenta. In that multi-channeling optimization procedure can be 

applied, which is incorporated in HAAG. In order to study the efficiency of this 
optimization we consider the same integration as before, but with all channels 
contributing to the generation and allowing for optimization. In this optimization 
procedure, we discard channels that contribute less than a certain pre-determined 
fraction to the set of available channels. It is expected, of course, that in end the 
right permutation will be 'chosen' by the optimization. This is indeed the case 
and the results are presented in Table 2. We see that the optimization results to 
a picture close to the one obtained with the single channel generation, with some 
noticeable improvement in the case of SARGE. 
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SARGE 
HAAG 


1 x 10 b 
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.483 x 10 _n 


.145 x 10~ 12 
.595 x 10" 13 


933 
18.2 


0.45 
1.19 



Table 2: Results for the all-channel generation with optimization. 



As is the case for any multi-channel generator, a computational complexity 
problem arises when the number of channels increases. For instance, in our case 
we are facing a number of \{n+ 1)! channels! On the other hand, it is also clear 
that the channels we are considering have a large overlap in most of the available 
phase space. It is therefore worth to investigate the dependence of the integration 
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efficiency on the number of channels used. This is presented in Table 3, where 
the full antenna 



S 2 Yl [ (Pl -P3)(P3 -P4)(P4 -P2)(P2 "PS) • • • (Pn+2 • Pi) ] 1 (24) 
P(2,...,n+2) 

is integrated, using a number of channels that has been selected on a random 
basis. We see the rather interesting phenomenon that a decent description can 
be achieved with a much smaller number of channels. Variations of this technique 



# channels 


2520 


1500 


1000 


500 


200 


50 


10 


/ 


5.33 


5.37 


5.48 


5.72 


6.14 


11.6 


84.7 


JVacc 
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26521 
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£(%) 


11.2 


13.1 


11.6 


7.1 


7.5 


1.7 


0.28 



Table 3: All-channel integration with subsets of channels for generation. 

of using only subsets of channels, for example choosing another subset after each 
step of multi-channel optimization, lead to the same picture. 

The complete results of the integration of the full antenna are presented in 
Table 4. We see that HAAG has a much better / factor than SARGE. On the other 
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SARGE 
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3 x 10 s 
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.176 x 10" ' 
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SARGE 
HAAG 


1.5 x 10 b 
1.8 x 10 5 


44529 
47911 


.157 x 10" ' 
.161 x 10" 7 


.135 x 10" M 
.905 x 10~ 10 


2.95 
7.15 
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5.7 
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SARGE 
HAAG 


1 x 10 7 
3.6 x 10 5 


47766 
45599 


.123 x lO" 7 
.123 x 10" 7 


.988 x 10- iu 
.241 x 10~ 10 


3.02 
5.11 


642 
13 
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SARGE 
HAAG 


1 x 10 8 
1 x 10 6 


53560 
49206 


.784 x 10" 8 
.789 x 10~ 8 


.554 x 10~ lu 
.496 x lO" 10 


3.29 
6.30 


4998 
39 



Table 4: Results for the all-channel integration. 



hand the e exhibits a less dramatic effect. This is related to the fact that SARGE 
generates a phase space that is much larger than the one defined by the cut on 
So- In that sense, if the main time consumption in a given computation is spent 
over the evaluation of the integrand (matrix element squared), it is more fair to 
compare the square of the estimated expected error, normalized by the number 
of accepted events iV acc . In that case we see that HAAG is still 2-3 times more 
efficient, and if we consider a smaller cut, namely ^/s^ =10 GeV, this gain goes 
up to an order of magnitude (Table 5). 

For a realistic calculation of the cross section of a QCD process, one may 
assume that the time it takes to perform one evaluation of the integrand is much 
larger than the time it takes to generate one accepted event and to calculate the 
weight. This means that the computing time is completely determined by the 
number of accepted events N acc . We introduce 
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Table 5: Results for the all-channel integrationwith s = 100 GeV 2 . 
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Figure 1: N^f /N gcn (a measure of computing time) as function of the number 
of produced partons. 



as a measure of the computing time. For a realistic calculation, one has to 
multiply this number by the evaluation time of the integrand, and divide by the 
square of the relative error one wants to reach. Fig.l shows this quantity as 
function of the number of produced partons using the data of Table 5. According 
to this graph, a calculation with SARGE would take 10 times longer than the 
calculation with HAAG. 

3 One-loop corrections to electroweak processes 

Scattering amplitudes in Quantum Field Theory can be represented by Feynman 
diagrams whose number grows extremely rapidly (faster than factorially) in the 
number of loops and external legs. This places severe limits in the calculation of 
multi-particle scattering amplitudes. The last few years, the development of new 
innovative methods and algorithms made possible to overcome these limitations. 
Based on Schwinger-Dyson recursive equations the complete scattering ampli- 
tude is computed directly without recourse to explicit Feynman diagrams. This 
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Y (4) Z(3) 



W-(1) Z(3) 



Z(3) T(4) 



f" "f 



Y(4) W+(2) 



Figure 2: Three of the six 1PI diagrams that contribute to the fermionic one- 
loop correction of the WWZ^ vertex. / represents a down-type fermion, and /' 
the corresponding up-type fermion. The other three diagrams are obtained by 
exchanging /«->/' and taking the opposite fermion current. 



results to a dramatic decrease in the computational cost which now depends only 
exponentially, i.e. ~ 3 n , in the number of external particles, n. Automatic com- 
putational tools, based on these Schwinger- Dyson recursive equations, have been 
developed that are able to describe any process within the Standard Electroweak 
theory and QCD [8, 9, 10, 11]. 

The high precision attainable by the future experiments calls for a reliable es- 
timate of the higher order corrections to the multi-particle scattering processes. 
This means that the full one-loop contributions along with the higher order lead- 
ing QED and electroweak corrections are necessary. Moreover, taking into ac- 
count the unstable particle contributions, re-summed propagator corrections have 
to be included. 

As a first step towards the extension of HELAC [8] towards the full one-loop 
level, re-summed boson propagators and fermion-loop corrections to boson ver- 
tices can be included. The main reasons to choose his collection of corrections 
is that it ensures gauge invariance [13], and that it is fairly straightforward to 
implement in an automatic program based on the Schwinger-Dyson method. If 
one, for example, wants to analyze processes that do not involve more than 6 
fermions, the one-loop corrections do not involve diagrams with more than 4 ex- 
ternal legs. So the "most complicated" diagrams that have to be included are 
1PI fermion-one-loop four-point functions (Fig.2). Below, we give a sketch how 
to evaluate them in a straightforward manner. The lower-point functions can be 
evaluated analogously. 

3.1 Evaluation of the one- loop four-point function 

Each diagram represents a function of the masses TOj of the fermions, the momenta 
Pi of the vector bosons, their polarization vectors Wi and their couplings Vi, ai to 
the fermion current. We conveniently choose to access them by permutating the 
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input of the general fermionic one-loop four-point function 

d n q 



= J -|Tr[P 1 (g)V 1 P 2 (g)V 2 P3(g)V 3 P4(g)V 4 ] , (26) 



with 



PM = , ^S = \t l + m 2 '. V, = M Vj + a jl5 ) . (27) 

[q + Ei=i Pi) - m j + 16 

Possible divergences are treated within the formalism of dimensional regulariza- 
tion. Ambiguities regarding 75 are avoided if the vectors pi and Wi are considered 
to be strictly 4-dimensional. The trace can be calculated with the help of com- 
puter algebra, for example with the program FORM [16], leading to 

f d n q num 4 (g; p li2 , 3 , "11,2,3,4, ^1,2,3,4, ^1,2,3,4, 01,2,3,4) , oc x 
J den 4 (g;pi, 2 ,3,m 1 , 2 ,3, 4 ) 

where the denominator function is defined by 



1 r. J-l , 2 



deni(q;p h2 ,...,i-i,m h 2,...,i) = J| + -m 2 j+i£ 

j=i 1 i=i 



(29) 



~~a — — ~ — v ' ( 30 ) 

m z dem 0; »i 9^,mi osa 



and the numerator num 4 is a fourth-order polynomial in the components of q. The 
integration problem is now reduced to that of the calculation of tensor integrals 
of the type 

d n q q vi --- q v 

iTT 2 den 4 (g;pi, 2 ,3,m 1 , 2 , 3 , 4 ) 

with r = 0,1,2,3,4. The Passarino-Veltman method [14, 15] uses the Lorentz 
covariance of these integrals to express them in terms of coefficient functions D 
through 

3 

A/i-i/ r (pi,2,3; "11,2,3,4) = Ph,vi • Pir.fr Aii 2 -i r (Pl,2,3; "11,2,3,4) , (31) 

ii,...,i r =0 

with D ili2 ... ir = if an odd number of indices are equal to 0, and the interpretation 

(32) 

The identification of the "D-functions and their expression in terms of the D- 
functions can also easily be performed by FORM. The output of FORM will contain 
symbols representing scalar products of the external momenta and the polar- 
ization vectors, contractions of these with the Levi-Civita symbol, and the IP- 
functions. This output can easily be turned into a FORTRAN-code. The scalar 
products and the Levi-Civita symbol are easy to be implemented, and the IP- 
functions can be extracted from the LoopTools-package [17], or one can make 
the effort to extend the FF-package [18], upon which LoopTools is based. 

2 We define the sum Y^d=i x i as a sum °f zer0 terms. 



11 



3.2 Result 



We will present the result of an actual calculation now. For more results and 
more details about the program presented above, we refer to [20]. There, we 
will also digress more about the well-known problem that the Passarino-Veltman 
method to calculate tensor integrals is numerically unstable for certain values of 
the external momenta. It involves the inversion of a kinematic matrix, which can 
become singular, although the one-loop function is perfectly-well defined. On 
might hope that, in a Monte Carlo calculation of a cross section, the probability 
to get too close to these phase-space points is too small to be concerned about 
this problem. We experienced, that this hope may be trusted for cross section 
calculations concerning processes at the coming generation of accelerators if the 
computation of the coefficient functions is performed at quadrupole precision. 
Restricting the use of quadrupole precision like this, the cpu-time stays within 
acceptable limits. We calculated the total cross section for the process 

e _ e + — > fj,~ v^u dr~ r + 

using the following cuts: Ei,E q > 5GeV for lepton and quark energies, a maximal 
cosine of 0.985 between all (initial and final state) charged leptons and quarks, 
and mu,m qq > lOGeV for the invariant masses of charged leptons and quarks. 
We used the renormalization scheme of [19], and the following input parameters 

m w = 80.35 GeV , m z = 91.1867 GeV 
Re[a (5) (7n^) _1 ] = 128.89 , a(0) _1 = 137.03599976 (33) 
G F = 1.16639 x 10~ 5 GeV . 

As far as the tree-order cross section is concerned we use the widely used Fixed 
Width scheme, where a fixed W^-boson width is implemented in all W^-boson 
propagators and where the G^-scheme is applied for evaluating the weak param- 
eters. We recall that the latter is defined by using m w , m z and G F as input 
parameters, together with the two relations 

2 ■> m w V2 2 2 

s w = 1 -—r > a = — G F m w s w ■ 

For the W and Z widths we use 

T w = 2.042GeV , T z = 2.49GeV . 

The results for E = 500GeV are a /ab = 54.96(26) a x jah = 57.31(28), and the 
K-factor is K/100 = 4.28(2). They show the expected contribution at the percent 
level of the higher order FL corrections to the total cross section. For comparison, 
the value of a /ab computed with 'Mega/WHIZARD [11, 12] is given by 55.07(19). 
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4 Conclusions 



We presented the algorithm HAAG, which uses the "hierarchical" strategy for phase 
space generation in order to efficiently map the antenna momentum structures, 
typically occurring in QCD amplitudes. It exhibits an improved efficiency com- 
pared to SARGE for multi-parton calculations, and it is more powerful in describing 
densities where a partial symmetrization over the permutation space is consid- 
ered. Also, HAAG makes no fundamental distinction among massless and massive 
particles, so it can be used for an arbitrary multi-partonic process. 

Furthermore we implemented one-loop corrections following the Fermion-Loop 
scheme in the program HELAC for automatic amplitude calculation and presented 
a result of a cross section calculation. 
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